C
      SUBROUTINE OPENFL(IN,ISTAT,FLNAME)
C ******************************************************************
C THIS SUBROUTINE OPENS AN INPUT/OUTPUT FILE ASSOCIATED WITH
C UNIT [IN], STATUS [ISTAT], AND FILE NAME [FLNAME].
C FILE IS OPENED AS 'FORMATTED' IF [IN]>0; 'UNFORMATTED' IF [IN]<0.
C FILE STATUS IS 'OLD' IF [ISTAT]>0; 'NEW' IF <0; 'UNKNOWN' IF =0.
C ******************************************************************
C last modified: 05-27-96
C

      IMPLICIT  NONE
      INTEGER   IN,ISTAT,I,IUNIT
      CHARACTER FLNAME*80,FLFORM*15,FLSTAT*15
	INTEGER*4 RESULT
C
C--DETERMINE FILE FORM AND STATUS
      IF(IN.GT.0) THEN
        FLFORM='FORMATTED'
      ELSEIF(IN.LT.0) THEN
c emrl        FLFORM='UNFORMATTED'
        FLFORM='BINARY'
      ELSE
        WRITE(*,100)
C-------EMRL JIG
        call stopfile
C-------EMRL JIG
        STOP
      ENDIF
      IF(ISTAT.GT.0) THEN
        FLSTAT='OLD'
      ELSEIF(ISTAT.LT.0) THEN
        FLSTAT='NEW'
      ELSEIF(ISTAT.EQ.0) THEN
        FLSTAT='UNKNOWN'
      ENDIF
  100 FORMAT(/1X,'ERROR: FILE CANNOT BE OPENED ON UNIT 0.')
C
C--GET FILE NAME IF NO DEFAULT NAME GIVEN
C     IF(IDFL.GT.0) GOTO 200
C
C     IF(FINDEX.EQ.' ') THEN
C       IF(IN.GT.0) WRITE(*,101) IN
C       IF(IN.LT.0) WRITE(*,102) -IN
C     ELSE
C       WRITE(*,103) FINDEX
C     ENDIF
C  10 READ(*,'(A50)') FLNAME
C     IF(FLNAME.EQ.' ') THEN
C       WRITE(*,*) 'Error: File Name Not Given.'
C       WRITE(*,*) 'Please Try Again =>'
C       GOTO 10
C     ENDIF
C 101 FORMAT(1X,'Enter Name of Formatted File for Unit',I3,': ')
C 102 FORMAT(1X,'Enter Name of Unformatted File for Unit',I3,': ')
C 103 FORMAT(1X,'Enter Name for ',A30)
C
C--OPEN FILE
  200 I=INDEX(FLNAME,' ')-1
      I=81
  201	I=I-1
      IF (FLNAME(I:I).EQ.' '.AND.I.GT.0) GOTO 201
      IUNIT = ABS(IN)
      OPEN(IUNIT,FILE=FLNAME(1:I),ERR=20,
     &  FORM=FLFORM,STATUS=FLSTAT)
      GOTO 30
C  20 IF(IDFL.GT.0) THEN
C       WRITE(*,2) FLNAME
C       STOP
C     ELSE
C       WRITE(*,3) FLNAME
C       GOTO 10
C     ENDIF
   20 WRITE(*,2) FLNAME
    2 FORMAT(1X,'Error: File Cannot Be Found or Opened =>',A80)
C   3 FORMAT(1X,'Error: File Cannot Be Found or Opened =>',A50
C    &      /1X,'Please Try Again =>')
C
C--FILE OPENED SUCCESSFUL, REWIND AND RETURN
   30 REWIND (IUNIT)
      RETURN
      END
C
C
      SUBROUTINE IARRAY(IA,ANAME,II,JJ,K,IN,IOUT)
C ************************************************************
C THIS SUBROUTINE IS USED TO INPUT 1 OR 2D INTEGER ARRAYS
C BY BLOCK, ZONAL, LIST-DIRECTED, UNFORMATTED,
C OR ANY USER-SPECIFIED FORMAT.
C ************************************************************
C last modified: 05-27-96
C

      IMPLICIT  NONE
      INTEGER   NZMAX
      PARAMETER (NZMAX=200)
      INTEGER   IA,II,JJ,K,IN,IOUT,IZV,IREAD,ICONST,IPRN,
     &          NBLOCK,NZONES,I1,I2,J1,J2,IZ,NN,I,J,N
      LOGICAL   OPD
      CHARACTER ANAME*24,FMTIN*20,FINDEX*30,FLNAME*50
      DIMENSION IA(JJ,II),IZV(NZMAX)
	INTEGER*4 RESULT
C
C--READ ARRAY CONTROL RECORD
C  =========================
      READ (IN,1) IREAD,ICONST,FMTIN,IPRN
    1 FORMAT(I10,I10,A20,I10)
C
C--IF IREAD=0, SET ALL ARRAY VALUES EQUAL TO ICONST.
C  =================================================
      IF(IREAD.NE.0) GOTO 50
C
      DO 10 I=1,II
        DO 12 J=1,JJ
          IA(J,I)=ICONST
   12   CONTINUE
   10 CONTINUE
      IF(K.GT.0) WRITE(IOUT,14) ANAME,ICONST,K
   14 FORMAT(39X,A24,' =',I15,' FOR LAYER',I3)
      IF(K.LE.0) WRITE(IOUT,16) ANAME,ICONST
   16 FORMAT(39X,A24,' =',I15)
      GOTO 500
C
C--IF IREAD=100, INPUT ARRAY USING FORMAT FMTIN
C  ============================================
   50 IF(IREAD.NE.100) GOTO 90
C
      IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IN,FMTIN
      IF(K.LE.0) WRITE(IOUT,22) ANAME,IN,FMTIN
   20 FORMAT(/21X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/21X,90('-'))
   22 FORMAT(/27X,A24,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/27X,77('-'))
      DO 30 I=1,II
        READ (IN,FMTIN) (IA(J,I),J=1,JJ)
   30 CONTINUE
      GOTO 300
C
C--IF IREAD=101, INPUT ARRAY USING BLOCK FORMAT
C  ============================================
   90 IF(IREAD.NE.101) GOTO 100
C
      IF(K.GT.0) WRITE(IOUT,55) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,60) ANAME,IN
   55 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/29X,72('-'))
   60 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/35X,59('-'))
C
C--READ NUMBER OF BLOCKS
      READ(IN,*) NBLOCK
C
C--READ VALUE OF EACH BLOCK
C--AND ASSIGN VALUE TO CELLS WITHIN THE BLOCK
      DO 70 N=1,NBLOCK
        READ(IN,*) I1,I2,J1,J2,IZ
        DO 72 I=I1,I2
          DO 74 J=J1,J2
            IA(J,I)=IZ
   74     CONTINUE
   72   CONTINUE
   70 CONTINUE
      GOTO 300
C
C--IF IREAD=102, INPUT ARRAY USING ZONAL FORMAT
C  ============================================
  100 IF(IREAD.NE.102) GOTO 200
C
      IF(K.GT.0) WRITE(IOUT,150) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,160) ANAME,IN
  150 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/29X,72('-'))
  160 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/35X,59('-'))
C
C--READ NUMBER OF ZONES
      READ(IN,*) NZONES
      IF(NZONES.GT.NZMAX) THEN
        WRITE(*,165)
C-------EMRL JIG
        call stopfile
C-------EMRL JIG
        STOP
      ENDIF
  165 FORMAT(/1X,'ERROR: MAXIMUM NUMBER OF ZONES EXCEEDED'
     & /1X,'INCREASE DIMENSION OF [NZMAX] IN SUBROUTINE [IARRAY]')
C
C--READ ZONAL MAP WITH FORMAT FMTIN
      READ(IN,*) (IZV(N),N=1,NZONES)
      DO 175 I=1,II
        READ(IN,FMTIN) (IA(J,I),J=1,JJ)
  175 CONTINUE
C
C--ASSIGN ZONAL VALUES
      DO 176 I=1,II
        DO 177 J=1,JJ
          NN=IA(J,I)
          IF(NN.EQ.0) THEN
            IA(J,I)=0
          ELSE
            IA(J,I)=IZV(NN)
          ENDIF
  177   CONTINUE
  176 CONTINUE
      GOTO 300
C
C--IF IREAD=103, INPUT ARRAY USING FREE FORMAT
C  ===========================================
  200 IF(IREAD.NE.103) GOTO 250
C
      IF(K.GT.0) WRITE(IOUT,210) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,220) ANAME,IN
  210 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/29X,71('-'))
  220 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/35X,58('-'))
C
C--READ ARRAY VALUES WITH FREE FORMAT
      READ(IN,*) ((IA(J,I),J=1,JJ),I=1,II)
      GOTO 300
C
C--IF IREAD IS EQUAL TO ANY OTHER VALUES,
C--READ ARRAY VALUES FROM AN EXTERNAL FILE ON UNIT [IREAD]
C  =======================================================
C
C--CHECK IF THE EXTERNAL FILE HAS BEEN OPENED.  IF NOT, OPEN
  250 INQUIRE(UNIT=IABS(IREAD),OPENED=OPD)
      IF(.NOT.OPD) THEN
        FINDEX=' '
        CALL OPENFL(IREAD,1,FLNAME)
      ENDIF
C
C--IF IREAD<0, INPUT ARRAY FROM AN UNFORMATTED FILE OM UNIT [-IREAD]
      IF(IREAD.LT.0) THEN
        IF(K.GT.0) WRITE(IOUT,256) ANAME,K,-IREAD
        IF(K.LE.0) WRITE(IOUT,258) ANAME,-IREAD
  256   FORMAT(/33X,A24,' FOR LAYER',I3,' READ UNFORMATTED',
     &  ' ON UNIT',I3/33X,65('-'))
  258   FORMAT(/40X,A24,' READ UNFORMATTED ON UNIT',
     &  I3/40X,52('-'))
C
C--READ AN UNFORMATTED DUMMY RECORD FIRST.
        READ(-IREAD)
        READ(-IREAD) IA
      ELSE
C
C--ELSE IF IREAD>0, INPUT ARRAY FROM AN FORMATTED FILE ON UNIT [IREAD]
C--WITH FORMAT FMTIN
        IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IREAD,FMTIN
        IF(K.LE.0) WRITE(IOUT,22) ANAME,IREAD,FMTIN
        DO 270 I=1,II
          READ (IREAD,FMTIN) (IA(J,I),J=1,JJ)
  270   CONTINUE
      ENDIF
C
C--IF ICONST NOT ZERO THEN MULTIPLY ARRAY VALUES BY ICONST.
  300 IF(ICONST.EQ.0) GO TO 320
      DO 302 I=1,II
      DO 304 J=1,JJ
        IA(J,I)=IA(J,I)*ICONST
  304 CONTINUE
  302 CONTINUE
C
C--IF PRINT CODE (IPRN) =>0 THEN PRINT ARRAY VALUES
  320 IF(IPRN.LT.0) GOTO 500
      CALL IPRINT(IA,ANAME,0,0,0,JJ,II,0,IPRN,IOUT)
C
C--RETURN
  500 RETURN
      END
C
C
      SUBROUTINE RARRAY(A,ANAME,II,JJ,K,IN,IOUT)
C ********************************************************
C THIS SUBROUTINE IS USED TO INPUT 1 OR 2D REAL ARRAYS,
C BY BLOCK, ZONAL, LIST-DIRECTED, UNFORMATTED,
C OR ANY USER-SPECIFIED FORMAT.
C ********************************************************
C last modified: 05-27-96
C

      IMPLICIT  NONE
      INTEGER   NZMAX
      PARAMETER (NZMAX=200)
      INTEGER   I,J,N,II,JJ,K,IN,IOUT,IREAD,IPRN,
     &          NBLOCK,NZONES,I1,I2,J1,J2,NN
      REAL      A,ZV,CONST,ZZ
      LOGICAL   OPD
      CHARACTER ANAME*24,FMTIN*20,FINDEX*30,FLNAME*50
      DIMENSION A(JJ,II),ZV(NZMAX)
	INTEGER*4 RESULT
C
C--READ ARRAY CONTROL RECORD
C  =========================
      READ (IN,1) IREAD,CONST,FMTIN,IPRN
    1 FORMAT(I10,F10.0,A20,I10)
C
C--IF IREAD=0, SET ALL ARRAY VALUES EQUAL TO CONST.
C  ================================================
      IF(IREAD.NE.0) GOTO 50
C
      DO 10 I=1,II
        DO 12 J=1,JJ
          A(J,I)=CONST
   12   CONTINUE
   10 CONTINUE
      IF(K.GT.0) WRITE(IOUT,14) ANAME,CONST,K
   14 FORMAT(39X,A24,' =',G15.7,' FOR LAYER',I3)
      IF(K.LE.0) WRITE(IOUT,16) ANAME,CONST
   16 FORMAT(39X,A24,' =',G15.7)
      GOTO 500
C
C--IF IREAD=100, INPUT ARRAY USING FORMAT FMTIN
C  ============================================
   50 IF(IREAD.NE.100) GOTO 90
C
      IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IN,FMTIN
      IF(K.LE.0) WRITE(IOUT,22) ANAME,IN,FMTIN
   20 FORMAT(/21X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/21X,90('-'))
   22 FORMAT(/27X,A24,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/27X,77('-'))
      DO 30 I=1,II
        READ (IN,FMTIN) (A(J,I),J=1,JJ)
   30 CONTINUE
      GOTO 300
C
C--IF IREAD=101, INPUT ARRAY USING BLOCK FORMAT
C  ============================================
   90 IF(IREAD.NE.101) GOTO 100
C
      IF(K.GT.0) WRITE(IOUT,55) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,60) ANAME,IN
   55 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/29X,72('-'))
   60 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/35X,59('-'))
C
C--READ NUMBER OF BLOCKS
      READ(IN,*) NBLOCK
C
C--READ VALUE OF EACH BLOCK
C--AND ASSIGN VALUE TO CELLS WITHIN THE BLOCK
      DO 70 N=1,NBLOCK
        READ(IN,*) I1,I2,J1,J2,ZZ
        DO 72 I=I1,I2
          DO 74 J=J1,J2
            A(J,I)=ZZ
   74     CONTINUE
   72   CONTINUE
   70 CONTINUE
      GOTO 300
C
C--IF IREAD=102, INPUT ARRAY USING ZONAL FORMAT
C  ============================================
  100 IF(IREAD.NE.102) GOTO 200
C
      IF(K.GT.0) WRITE(IOUT,150) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,160) ANAME,IN
  150 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/29X,72('-'))
  160 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/35X,59('-'))
C
C--READ NUMBER OF ZONES
      READ(IN,*) NZONES
      IF(NZONES.GT.NZMAX) THEN
        WRITE(*,165)
C-------EMRL JIG
        call stopfile
C-------EMRL JIG
        STOP
      ENDIF
  165 FORMAT(1X,'ERROR: MAXIMUM NUMBER OF ZONES EXCEEDED'
     & /1X,'INCREASE DIMENSION OF [NZMAX] IN SUBROUTINE [RARRAY]')
C
C--READ ZONAL MAP WITH FORMAT FMTIN
      READ(IN,*) (ZV(N),N=1,NZONES)
      DO 175 I=1,II
        READ(IN,FMTIN) (A(J,I),J=1,JJ)
  175 CONTINUE
C
C--ASSIGN ZONAL VALUES
      DO 176 I=1,II
        DO 177 J=1,JJ
          NN=A(J,I)
          IF(NN.EQ.0) THEN
            A(J,I)=0
          ELSE
            A(J,I)=ZV(NN)
          ENDIF
  177   CONTINUE
  176 CONTINUE
      GOTO 300
C
C--IF IREAD=103, INPUT ARRAY USING FREE FORMAT
C  ===========================================
  200 IF(IREAD.NE.103) GOTO 250
C
      IF(K.GT.0) WRITE(IOUT,210) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,220) ANAME,IN
  210 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/29X,71('-'))
  220 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/35X,58('-'))
C
C--READ ARRAY VALUES WITH FREE FORMAT
      READ(IN,*) ((A(J,I),J=1,JJ),I=1,II)
      GOTO 300
C
C--IF IREAD IS EQUAL TO ANY OTHER VALUES,
C--READ ARRAY VALUES FROM AN EXTERNAL FILE ON UNIT [IREAD]
C  =======================================================
C
C--CHECK IF THE EXTERNAL FILE HAS BEEN OPENED.  IF NOT, OPEN
  250 INQUIRE(UNIT=IABS(IREAD),OPENED=OPD)
      IF(.NOT.OPD) THEN
        FINDEX=' '
        CALL OPENFL(IREAD,1,FLNAME)
      ENDIF
C
C--IF IREAD<0, INPUT ARRAY FROM AN UNFORMATTED FILE OM UNIT [-IREAD]
      IF(IREAD.LT.0) THEN
        IF(K.GT.0) WRITE(IOUT,256) ANAME,K,-IREAD
        IF(K.LE.0) WRITE(IOUT,258) ANAME,-IREAD
  256   FORMAT(/33X,A24,' FOR LAYER',I3,' READ UNFORMATTED',
     &  ' ON UNIT',I3/33X,65('-'))
  258   FORMAT(/40X,A24,' READ UNFORMATTED ON UNIT',
     &  I3/40X,52('-'))
C
C--READ AN UNFORMATTED DUMMY RECORD FIRST.
        READ(-IREAD)
        READ(-IREAD) A
      ELSE
C
C--ELSE IF IREAD>0, INPUT ARRAY FROM AN FORMATTED FILE ON UNIT [IREAD]
C--WITH FORMAT FMTIN
        IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IREAD,FMTIN
        IF(K.LE.0) WRITE(IOUT,22) ANAME,IREAD,FMTIN
        DO 270 I=1,II
          READ (IREAD,FMTIN) (A(J,I),J=1,JJ)
  270   CONTINUE
      ENDIF
C
C--IF CONST NOT ZERO THEN MULTIPLY ARRAY VALUES BY CONST.
  300 IF(CONST.EQ.0) GO TO 320
      DO 302 I=1,II
      DO 304 J=1,JJ
        A(J,I)=A(J,I)*CONST
  304 CONTINUE
  302 CONTINUE
C
C--IF PRINT CODE (IPRN) =>0 THEN PRINT ARRAY VALUES.
  320 IF(IPRN.LT.0) RETURN
      CALL RPRINT(A,ANAME,0,0,0,JJ,II,0,IPRN,IOUT)
C
C8------RETURN
  500 RETURN
      END
C
C
      SUBROUTINE IPRINT(IA,TEXT,KTRN,KSTP,KPER,NCOL,NROW,
     & ILAY,IPRN,IOUT)
C ************************************************************
C PRINT AN INTEGER 1 OR 2D ARRAY IN WRAP OR STRIP FORM.
C [MODIFIED FROM MCDONALD AND HARBAUGH (1988)].
C ************************************************************
C last modified: 05-27-96
C
      IMPLICIT  NONE
      INTEGER   IA,KTRN,KSTP,KPER,NCOL,NROW,ILAY,IPRN,IOUT,
     &          IP,NCPF,NCAP,J1,J2,NSTRIP,ISP,J,I,N
      CHARACTER TEXT*16
      DIMENSION IA(NCOL,NROW)
C
C--PRINT A HEADER
      IF(ILAY.LE.0) GO TO 5
      IF(KTRN.GT.0) GO TO 4
      WRITE(IOUT,1) TEXT,ILAY,KSTP,KPER
    1 FORMAT(/35X,A16,' IN LAYER',I3,
     & ' FOR TIME STEP',I3,', STRESS PERIOD',I3/35X,63('-'))
      GOTO 5
C
    4 WRITE(IOUT,2) TEXT,ILAY,KTRN,KSTP,KPER
    2 FORMAT(/21X,A16,' IN LAYER',I3,' AT END OF TRANSPORT STEP',I5,
     & ', TIME STEP',I3,', STRESS PERIOD',I3/21X,90('-'))
C
C--MAKE SURE IPRN VALUE IS WITHIN PRINT FORMAT-CODE RANGE
    5 IP=IPRN
      IF(IP.GT.5.OR.IP.LT.-5) IP=0
      IF(IP.GE.0) IP=IP+1
      IF(IP.LT.0) IP=IP-1
C
C--DETERMINE THE NUMBER OF VALUES (NCAP) PRINTED ON ONE LINE,
C--NUMBER OF BLANK SPACES TO LEAVE AT START OF THE LINE (ISP)
C--AND NUMBER OF SPACES IN EACH COLUMN FIELD (NCPF)
      IF(IABS(IP).EQ.1) THEN
        NCPF=12
        ISP=4
        NCAP=10
      ELSEIF(IABS(IP).GT.1) THEN
        NCPF=IABS(IP)
        ISP=4
        NCAP=125/IABS(IP)/5*5
      ENDIF
C
C--IF IP>0, ARRAY IS PRINTED IN WRAP FORMAT
      IF(IP.GT.0.OR.NCOL.LE.NCAP) THEN
        NSTRIP=1
        J1=1
        J2=NCOL
C
C--ELSE IF IP<0, ARRAY IS PRINTED IN STRIP FORMAT.
C--NUMBER OF STRIPS IS CALCULATED AS [NSTRIP]
      ELSE
        NSTRIP=(NCOL-1)/NCAP+1
        J1=1-NCAP
        J2=0
      ENDIF
C
C--LOOP THROUGH THE STRIPS
      DO 400 N=1,NSTRIP
C
C--CALCULATE FIRST(J1) & LAST(J2) COLUMNS FOR THIS STRIP
C--IF STRIP FORM IS USED
      IF(NSTRIP.GT.1) THEN
        J1=J1+NCAP
        J2=J2+NCAP
        IF(J2.GT.NCOL) J2=NCOL
      ENDIF
C
C--PRINT COLUMN NUMBERS ABOVE THE STRIP.
      CALL UCOLNO(J1,J2,ISP,NCAP,NCPF,IOUT)
C
C--PRINT EACH ROW
      DO 410 I=1,NROW
C
C--SELECT THE FORMAT
        GOTO (401,402,403,404,405,406), IABS(IP)
C
C--FORMAT 10I11
  401   IF(IP.GT.0) WRITE(IOUT,1001) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2001) I,(IA(J,I),J=J1,J2)
 1001   FORMAT(1X,I3,2X,I11,9(1X,I11)/(5X,10(1X,I11)))
 2001   FORMAT(1X,I3,2X,I11,9(1X,I11))
        GO TO 410
C
C--FORMAT 60I1
  402   IF(IP.GT.0) WRITE(IOUT,1002) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2002) I,(IA(J,I),J=J1,J2)
 1002   FORMAT(1X,I3,1X,60(1X,I1)/(5X,60(1X,I1)))
 2002   FORMAT(1X,I3,1X,60(1X,I1))
        GO TO 410
C
C--FORMAT 40I2
  403   IF(IP.GT.0) WRITE(IOUT,1003) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2003) I,(IA(J,I),J=J1,J2)
 1003   FORMAT(1X,I3,1X,40(1X,I2)/(5X,40(1X,I2)))
 2003   FORMAT(1X,I3,1X,40(1X,I2))
        GO TO 410
C
C--FORMAT 30I3
  404   IF(IP.GT.0) WRITE(IOUT,1004) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2004) I,(IA(J,I),J=J1,J2)
 1004   FORMAT(1X,I3,1X,30(1X,I3)/(5X,30(1X,I3)))
 2004   FORMAT(1X,I3,1X,30(1X,I3))
        GO TO 410
C
C--FORMAT 25I4
  405   IF(IP.GT.0) WRITE(IOUT,1005) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2005) I,(IA(J,I),J=J1,J2)
 1005   FORMAT(1X,I3,1X,25(1X,I4)/(5X,25(1X,I4)))
 2005   FORMAT(1X,I3,1X,25(1X,I4))
        GO TO 410
C
C--FORMAT 20I5
  406   IF(IP.GT.0) WRITE(IOUT,1006) I,(IA(J,I),J=J1,J2)
        IF(IP.LT.0) WRITE(IOUT,2006) I,(IA(J,I),J=J1,J2)
 1006   FORMAT(1X,I3,1X,20(1X,I5)/(5X,20(1X,I5)))
 2006   FORMAT(1X,I3,1X,20(1X,I5))
  410 CONTINUE
C
  400 CONTINUE
C
C--RETURN
      RETURN
      END
C
C
      SUBROUTINE RPRINT(BUFF,TEXT,KTRN,KSTP,KPER,NCOL,NROW,
     & ILAY,IPRN,IOUT)
C ****************************************************************
C PRINT A REAL 1 OR 2D ARRAY IN WRAP OR STRIP FORM.
C [MODIFIED FROM MCDONALD AND HARBAUGH (1988)].
C ****************************************************************
C last modified: 05-27-96
C
      IMPLICIT  NONE
      INTEGER   KTRN,KSTP,KPER,NCOL,NROW,ILAY,IPRN,IOUT,IP,J,I,
     &          J1,J2,NCPF,NCAP,NSTRIP,ISP,N
      REAL      BUFF
      CHARACTER TEXT*16
      DIMENSION BUFF(NCOL,NROW)
C
C--PRINT A HEADER
      IF(ILAY.LE.0) GO TO 5
      IF(KTRN.GT.0) GO TO 4
      WRITE(IOUT,1) TEXT,ILAY,KSTP,KPER
    1 FORMAT(/35X,A16,' IN LAYER',I3,
     & ' FOR TIME STEP',I3,', STRESS PERIOD',I3/35X,63('-'))
      GOTO 5
C
    4 WRITE(IOUT,2) TEXT,ILAY,KTRN,KSTP,KPER
    2 FORMAT(/21X,A16,' IN LAYER',I3,' AT END OF TRANSPORT STEP',I5,
     & ', TIME STEP',I3,', STRESS PERIOD',I3/21X,90('-'))
C
C--MAKE SURE IPRN VALVE IS WITHIN PRINT-FORMAT CODE RANGE
    5 IP=IPRN
      IF(IP.GT.12) IP=12
      IF(IP.LT.-12) IP=-12
      IF(IP.EQ.0) IP=12
C
C--DETERMINE THE NUMBER OF VALUES (NCAP) PRINTED ON ONE LINE,
C--NUMBER OF BLANK SPACES TO LEAVE AT START OF THE LINE (ISP)
C--AND NUMBER OF SPACES IN EACH COLUMN FIELD (NCPF)
      IF(IABS(IP).EQ.1) NCAP=11
      IF(IABS(IP).EQ.2) NCAP=9
      IF(IABS(IP).GT.2 .AND. IABS(IP).LT.7) NCAP=15
      IF(IABS(IP).GT.6 .AND. IABS(IP).LT.12) NCAP=20
      IF(IABS(IP).EQ.12) NCAP=10
      NCPF=129/NCAP
      ISP=0
      IF(NCAP.GT.12) ISP=3
C
C--IF IP>0, ARRAY IS PRINTED IN WRAP FORM
      IF(IP.GT.0.OR.NCOL.LE.NCAP) THEN
        NSTRIP=1
        J1=1
        J2=NCOL
C
C--ELSE IF IP<0, ARRAY IS PRINTED IN STRIP FORM.
C--NUMBER OF STRIPS IS CALCULATED AS [NSTRIP]
      ELSE
        NSTRIP=(NCOL-1)/NCAP + 1
        J1=1-NCAP
        J2=0
      ENDIF
C
C--LOOP THROUGH THE STRIPS
      DO 2000 N=1,NSTRIP
C
C--CALCULATE FIRST(J1) & LAST(J2) COLUMNS FOR THIS STRIP
C--IF STRIP FORM IS USED
      IF(NSTRIP.GT.1) THEN
        J1=J1+NCAP
        J2=J2+NCAP
        IF(J2.GT.NCOL) J2=NCOL
      ENDIF
C
C--PRINT COLUMN NUMBERS ABOVE THE STRIP
      CALL UCOLNO(J1,J2,ISP,NCAP,NCPF,IOUT)
C
C--LOOP THROUGH THE ROWS PRINTING COLS J1 THRU J2 WITH FORMAT IP
      DO 1000 I=1,NROW
      GO TO(10,20,30,40,50,60,70,80,90,100,110,120), IABS(IP)
C
C--FORMAT 11G10.3
   10 IF(IP.GT.0) WRITE(IOUT,11) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,12) I,(BUFF(J,I),J=J1,J2)
   11 FORMAT(1X,I3,2X,1PG10.3,10(1X,G10.3)/(5X,11(1X,G10.3)))
   12 FORMAT(1X,I3,2X,1PG10.3,10(1X,G10.3))
      GO TO 1000
C
C--FORMAT 9G13.6
   20 IF(IP.GT.0) WRITE(IOUT,21) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,22) I,(BUFF(J,I),J=J1,J2)
   21 FORMAT(1X,I3,2X,1PG13.6,8(1X,G13.6)/(5X,9(1X,G13.6)))
   22 FORMAT(1X,I3,2X,1PG13.6,8(1X,G13.6))
      GO TO 1000
C
C--FORMAT 15F7.1
   30 IF(IP.GT.0) WRITE(IOUT,31) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,32) I,(BUFF(J,I),J=J1,J2)
   31 FORMAT(1X,I3,1X,15(1X,F7.1)/(5X,15(1X,F7.1)))
   32 FORMAT(1X,I3,1X,15(1X,F7.1))
      GO TO 1000
C
C--FORMAT 15F7.2
   40 IF(IP.GT.0) WRITE(IOUT,41) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,42) I,(BUFF(J,I),J=J1,J2)
   41 FORMAT(1X,I3,1X,15(1X,F7.2)/(5X,15(1X,F7.2)))
   42 FORMAT(1X,I3,1X,15(1X,F7.2))
      GO TO 1000
C
C--FORMAT 15F7.3
   50 IF(IP.GT.0) WRITE(IOUT,51) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,52) I,(BUFF(J,I),J=J1,J2)
   51 FORMAT(1X,I3,1X,15(1X,F7.3)/(5X,15(1X,F7.3)))
   52 FORMAT(1X,I3,1X,15(1X,F7.3))
      GO TO 1000
C
C--FORMAT 15F7.4
   60 IF(IP.GT.0) WRITE(IOUT,61) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,62) I,(BUFF(J,I),J=J1,J2)
   61 FORMAT(1X,I3,1X,15(1X,F7.4)/(5X,15(1X,F7.4)))
   62 FORMAT(1X,I3,1X,15(1X,F7.4))
      GO TO 1000
C
C--FORMAT 20F5.0
   70 IF(IP.GT.0) WRITE(IOUT,71) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,72) I,(BUFF(J,I),J=J1,J2)
   71 FORMAT(1X,I3,1X,20(1X,F5.0)/(5X,20(1X,F5.0)))
   72 FORMAT(1X,I3,1X,20(1X,F5.0))
      GO TO 1000
C
C--FORMAT 20F5.1
   80 IF(IP.GT.0) WRITE(IOUT,81) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,82) I,(BUFF(J,I),J=J1,J2)
   81 FORMAT(1X,I3,1X,20(1X,F5.1)/(5X,20(1X,F5.1)))
   82 FORMAT(1X,I3,1X,20(1X,F5.1))
      GO TO 1000
C
C--FORMAT 20F5.2
   90 IF(IP.GT.0) WRITE(IOUT,91) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,92) I,(BUFF(J,I),J=J1,J2)
   91 FORMAT(1X,I3,1X,20(1X,F5.2)/(5X,20(1X,F5.2)))
   92 FORMAT(1X,I3,1X,20(1X,F5.2))
      GO TO 1000
C
C--FORMAT 20F5.3
  100 IF(IP.GT.0) WRITE(IOUT,101) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,102) I,(BUFF(J,I),J=J1,J2)
  101 FORMAT(1X,I3,1X,20(1X,F5.3)/(5X,20(1X,F5.3)))
  102 FORMAT(1X,I3,1X,20(1X,F5.3))
      GO TO 1000
C
C--FORMAT 20F5.4
  110 IF(IP.GT.0) WRITE(IOUT,111) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,112) I,(BUFF(J,I),J=J1,J2)
  111 FORMAT(1X,I3,1X,20(1X,F5.4)/(5X,20(1X,F5.4)))
  112 FORMAT(1X,I3,1X,20(1X,F5.4)/(5X,20(1X,F5.4)))
      GO TO 1000
C
C--FORMAT 10G11.4
  120 IF(IP.GT.0) WRITE(IOUT,121) I,(BUFF(J,I),J=J1,J2)
      IF(IP.LT.0) WRITE(IOUT,122) I,(BUFF(J,I),J=J1,J2)
  121 FORMAT(1X,I3,2X,1PG11.4,9(1X,G11.4)/(5X,10(1X,G11.4)))
  122 FORMAT(1X,I3,2X,1PG11.4,9(1X,G11.4))
C
 1000 CONTINUE
 2000 CONTINUE
C
C--RETURN
      RETURN
      END
C
C
      SUBROUTINE UCOLNO(NLBL1,NLBL2,NSPACE,NCPL,NDIG,IOUT)
C ****************************************************************
C OUTPUT COLUMN NUMBERS ABOVE A MATRIX PRINTOUT.
C [MODIFIED FROM MCDONALD AND HARBAUGH (1988)].
C ****************************************************************
C last modified: 05-27-96
C
      IMPLICIT  NONE
      INTEGER   NLBL1,NLBL2,NSPACE,NCPL,NDIG,IOUT,N,
     &          J1,J2,J,I1,I2,I3,I,NTOT,NWRAP,NBF,NLBL
      CHARACTER DOT*4,SPACE*4,DG*4,BF*4
      DIMENSION BF(130),DG(10)
C
C--ASSIGN CHARACTER STRING
      DG(1)='0   '
      DG(2)='1   '
      DG(3)='2   '
      DG(4)='3   '
      DG(5)='4   '
      DG(6)='5   '
      DG(7)='6   '
      DG(8)='7   '
      DG(9)='8   '
      DG(10)='9   '
      DOT='.  '
      SPACE='    '
C
C--CALCULATE # OF COLUMNS TO BE PRINTED (NLBL), WIDTH
C--OF A LINE (NTOT), NUMBER OF LINES (NWRAP).
      WRITE(IOUT,1)
    1 FORMAT(1X)
      NLBL=NLBL2-NLBL1+1
      N=NLBL
      IF(NLBL.GT.NCPL) N=NCPL
      NTOT=NSPACE+N*NDIG
      IF(NTOT.GT.130) GO TO 50
      NWRAP=(NLBL-1)/NCPL + 1
      J1=NLBL1-NCPL
      J2=NLBL1-1
C
C--BUILD AND PRINT EACH LINE
      DO 40 N=1,NWRAP
C
C--CLEAR THE BUFFER (BF).
      DO 20 I=1,130
      BF(I)=SPACE
   20 CONTINUE
      NBF=NSPACE
C
C--DETERMINE FIRST (J1) AND LAST (J2) COLUMN # FOR THIS LINE.
      J1=J1+NCPL
      J2=J2+NCPL
      IF(J2.GT.NLBL2) J2=NLBL2
C
C--LOAD THE COLUMN #'S INTO THE BUFFER.
      DO 30 J=J1,J2
      NBF=NBF+NDIG
      I2=J/10
      I1=J-I2*10+1
      BF(NBF)=DG(I1)
      IF(I2.EQ.0) GO TO 30
      I3=I2/10
      I2=I2-I3*10+1
      BF(NBF-1)=DG(I2)
      IF(I3.EQ.0) GO TO 30
      BF(NBF-2)=DG(I3+1)
   30 CONTINUE
C
C--PRINT THE CONTENTS OF THE BUFFER (I.E. PRINT THE LINE).
      WRITE(IOUT,31) (BF(I),I=1,NBF)
   31 FORMAT(1X,130A1)
C
   40 CONTINUE
C
C--PRINT A LINE OF DOTS (FOR ESTHETIC PURPOSES ONLY).
   50 NTOT=NTOT+5
      IF(NTOT.GT.130) NTOT=130
      WRITE(IOUT,51) (DOT,I=1,NTOT)
   51 FORMAT(1X,130A1)
C
C--RETURN
      RETURN
      END

C
      SUBROUTINE DPRARRAY(A,ANAME,II,JJ,K,IN,IOUT)
C ********************************************************
C THIS SUBROUTINE IS USED TO INPUT 1 or 2D DOUBLE PRECISION ARRAYS,
C BY BLOCK, ZONAL, LIST-DIRECTED, UNFORMATTED,
C OR ANY USER-SPECIFIED FORMAT.
C ********************************************************
C last modified: 14-03-97
C

      IMPLICIT  NONE
      INTEGER   NZMAX
      PARAMETER (NZMAX=200)
      INTEGER   I,J,N,II,JJ,K,IN,IOUT,IREAD,IPRN,
     &          NBLOCK,NZONES,I1,I2,J1,J2,NN
      double precision  A,ZV,CONST,ZZ
      LOGICAL   OPD
      CHARACTER ANAME*24,FMTIN*20,FINDEX*30,FLNAME*50
      DIMENSION A(JJ,II),ZV(NZMAX)
	INTEGER*4 RESULT
C
C--READ ARRAY CONTROL RECORD
C  =========================
      READ (IN,1) IREAD,CONST,FMTIN,IPRN
    1 FORMAT(I10,F10.0,A20,I10)
C
C--IF IREAD=0, SET ALL ARRAY VALUES EQUAL TO CONST.
C  ================================================
      IF(IREAD.NE.0) GOTO 50
C
      DO 10 I=1,II
        DO 12 J=1,JJ
          A(J,I)=CONST
   12   CONTINUE
   10 CONTINUE
      IF(K.GT.0) WRITE(IOUT,14) ANAME,CONST,K
   14 FORMAT(39X,A24,' =',G15.7,' FOR LAYER',I3)
      IF(K.LE.0) WRITE(IOUT,16) ANAME,CONST
   16 FORMAT(39X,A24,' =',G15.7)
      GOTO 500
C
C--IF IREAD=100, INPUT ARRAY USING FORMAT FMTIN
C  ============================================
   50 IF(IREAD.NE.100) GOTO 90
C
      IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IN,FMTIN
      IF(K.LE.0) WRITE(IOUT,22) ANAME,IN,FMTIN
   20 FORMAT(/21X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/21X,90('-'))
   22 FORMAT(/27X,A24,' READ ON UNIT',
     & I3,' USING FORMAT: "',A20,'"'/27X,77('-'))
      DO 30 I=1,II
        READ (IN,FMTIN) (A(J,I),J=1,JJ)
   30 CONTINUE
      GOTO 300
C
C--IF IREAD=101, INPUT ARRAY USING BLOCK FORMAT
C  ============================================
   90 IF(IREAD.NE.101) GOTO 100
C
      IF(K.GT.0) WRITE(IOUT,55) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,60) ANAME,IN
   55 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/29X,72('-'))
   60 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING BLOCK FORMAT'/35X,59('-'))
C
C--READ NUMBER OF BLOCKS
      READ(IN,*) NBLOCK
C
C--READ VALUE OF EACH BLOCK
C--AND ASSIGN VALUE TO CELLS WITHIN THE BLOCK
      DO 70 N=1,NBLOCK
        READ(IN,*) I1,I2,J1,J2,ZZ
        DO 72 I=I1,I2
          DO 74 J=J1,J2
            A(J,I)=ZZ
   74     CONTINUE
   72   CONTINUE
   70 CONTINUE
      GOTO 300
C
C--IF IREAD=102, INPUT ARRAY USING ZONAL FORMAT
C  ============================================
  100 IF(IREAD.NE.102) GOTO 200
C
      IF(K.GT.0) WRITE(IOUT,150) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,160) ANAME,IN
  150 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/29X,72('-'))
  160 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING ZONAL FORMAT'/35X,59('-'))
C
C--READ NUMBER OF ZONES
      READ(IN,*) NZONES
      IF(NZONES.GT.NZMAX) THEN
        WRITE(*,165)
C-------EMRL JIG
        call stopfile
C-------EMRL JIG
        STOP
      ENDIF
  165 FORMAT(1X,'ERROR: MAXIMUM NUMBER OF ZONES EXCEEDED'
     & /1X,'INCREASE DIMENSION OF [NZMAX] IN SUBROUTINE [RARRAY]')
C
C--READ ZONAL MAP WITH FORMAT FMTIN
      READ(IN,*) (ZV(N),N=1,NZONES)
      DO 175 I=1,II
        READ(IN,FMTIN) (A(J,I),J=1,JJ)
  175 CONTINUE
C
C--ASSIGN ZONAL VALUES
      DO 176 I=1,II
        DO 177 J=1,JJ
          NN=A(J,I)
          IF(NN.EQ.0) THEN
            A(J,I)=0
          ELSE
            A(J,I)=ZV(NN)
          ENDIF
  177   CONTINUE
  176 CONTINUE
      GOTO 300
C
C--IF IREAD=103, INPUT ARRAY USING FREE FORMAT
C  ===========================================
  200 IF(IREAD.NE.103) GOTO 250
C
      IF(K.GT.0) WRITE(IOUT,210) ANAME,K,IN
      IF(K.LE.0) WRITE(IOUT,220) ANAME,IN
  210 FORMAT(/29X,A24,' FOR LAYER',I3,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/29X,71('-'))
  220 FORMAT(/35X,A24,' READ ON UNIT',
     & I3,' USING FREE FORMAT'/35X,58('-'))
C
C--READ ARRAY VALUES WITH FREE FORMAT
      READ(IN,*) ((A(J,I),J=1,JJ),I=1,II)
      GOTO 300
C
C--IF IREAD IS EQUAL TO ANY OTHER VALUES,
C--READ ARRAY VALUES FROM AN EXTERNAL FILE ON UNIT [IREAD]
C  =======================================================
C
C--CHECK IF THE EXTERNAL FILE HAS BEEN OPENED.  IF NOT, OPEN
  250 INQUIRE(UNIT=IABS(IREAD),OPENED=OPD)
      IF(.NOT.OPD) THEN
        FINDEX=' '
        CALL OPENFL(IREAD,1,FLNAME)
      ENDIF
C
C--IF IREAD<0, INPUT ARRAY FROM AN UNFORMATTED FILE OM UNIT [-IREAD]
      IF(IREAD.LT.0) THEN
        IF(K.GT.0) WRITE(IOUT,256) ANAME,K,-IREAD
        IF(K.LE.0) WRITE(IOUT,258) ANAME,-IREAD
  256   FORMAT(/33X,A24,' FOR LAYER',I3,' READ UNFORMATTED',
     &  ' ON UNIT',I3/33X,65('-'))
  258   FORMAT(/40X,A24,' READ UNFORMATTED ON UNIT',
     &  I3/40X,52('-'))
C
C--READ AN UNFORMATTED DUMMY RECORD FIRST.
        READ(-IREAD)
        READ(-IREAD) A
      ELSE
C
C--ELSE IF IREAD>0, INPUT ARRAY FROM AN FORMATTED FILE ON UNIT [IREAD]
C--WITH FORMAT FMTIN
        IF(K.GT.0) WRITE(IOUT,20) ANAME,K,IREAD,FMTIN
        IF(K.LE.0) WRITE(IOUT,22) ANAME,IREAD,FMTIN
        DO 270 I=1,II
          READ (IREAD,FMTIN) (A(J,I),J=1,JJ)
  270   CONTINUE
      ENDIF
C
C--IF CONST NOT ZERO THEN MULTIPLY ARRAY VALUES BY CONST.
  300 IF(CONST.EQ.0) GO TO 320
      DO 302 I=1,II
      DO 304 J=1,JJ
        A(J,I)=A(J,I)*CONST
  304 CONTINUE
  302 CONTINUE
C
C--IF PRINT CODE (IPRN) =>0 THEN PRINT ARRAY VALUES.
  320 IF(IPRN.LT.0) RETURN
C      CALL RPRINT(A,ANAME,0,0,0,JJ,II,0,IPRN,IOUT) !print commented out
C
C8------RETURN
  500 RETURN
      END

